X <- matrix(c(8,1,0,4,6,5,6,8,7,10,4,7,8,2,5,0,3,6), nrow = 6, ncol = 3, byrow
            = TRUE)

#----- Estudiantes--------
estudiantes_datos <- read.table('EjemploEstudiantes.csv', header=TRUE, sep=';',dec=',',row.names=1)
estudiantes_datos<- as.matrix(estudiantes_datos)

#----- beans--------------
beans_datos <- read.csv("beansV2.csv")

#Convertir Class a código disyuntivo
categorias_class  <- unique(beans_datos$Class)
lista_disyuntivo  <- list()
beans_datos_disyuntivo <- beans_datos[,-17]

for (i in 1:length(categorias_class)) {
  lista_disyuntivo[[i]] <- as.numeric(beans_datos$Class == categorias_class[i])
  names(lista_disyuntivo)[i] <- paste("Class_", categorias_class[i], sep = "")
  beans_datos_disyuntivo <- cbind(beans_datos_disyuntivo, lista_disyuntivo[i])
}
beans_datos <- as.matrix(beans_datos_disyuntivo)

#Centramos y reducimos la matriz con la siguiente función
centrar_y_reducir <-  function(matriz){
  # Se obtiene la cantidad de filas y columnas que posee la matriz.
  columnas <- ncol(matriz)
  filas <- nrow(matriz)
  
  # Se crea una nueva matriz con los datos centrados y reducidos.
  for(i in 1:columnas){
    matriz[,i] <- (matriz[,i] - mean(matriz[,i])) / (sd(matriz[,i])*sqrt((filas - 1)/filas))
  }
  
  return(matriz)
}

#Calcular la matriz de correlaciones 
matriz_correlaciones_ind <- function(matriz) {
  n <- nrow(matriz) 
  resultado <- (1/n)*t(matriz)%*%matriz
  return(resultado)
}

#Calcular la matriz de calidades de los individuos (cosenos cuadrados)
matriz_cos2_ind <- function(C, matriz) {
  n <-nrow(matriz)
  m <-ncol(matriz)
  resultado <- matrix(0, n,m)
  
  for(i in 1: n){
    suma <- 0 
    for(j in 1: m){
      suma <- suma + matriz[i,j]^2
    }
    for(r in 1: m){
      resultado[i,r] <- (C[i,r]^2)/suma
    }
  }
  return(resultado)
}

#Calcular la contribución de cada individuo i a la varianza total del eje r
contribucion_ind <- function(C, valores_propios) {
  n <-nrow(C)
  m <-ncol(C)
  resultado <- matrix(0, n, m)
  
  for(i in 1: n){
    for(r in 1: m){
      resultado[i,r] <- ((C[i,r]^2)/(n*valores_propios[r]))*100
    }
  }
  return(resultado)
}

# Calcular la matriz de coordenadas de las variables
coordenadas_var <- function(V, valores_propios){
  m <- ncol(V)
  resultado <- matrix(0, m, m)
  
  for(r in 1:m){
      resultado[,r] <- sqrt(valores_propios[r])*V[,r]
    }
    
  return(resultado)
}

#Obtener contribuciones de las variables al cálculo de la varianza total del eje r
contribuciones_var <- function(cos2){
  m <- ncol(cos2)
  contribuciones <- matrix(nrow = m, ncol = m)
  
  for(i in 1:m){
    for(j in 1:m){
      contribuciones[j,i] <- (cos2[j,i]/sum(cos2[,i]))*100
    }
  }
  
  return(contribuciones)
}

#Calcular vector de inercias de los ejes
calcular_inercias <- function(valores_propios){
  m <- length(valores_propios)
  resultado <- c()
  
  for(j in 1: m){
    resultado[j] <-100*(valores_propios[j]/m)
  }
  return(resultado)
}

ACP_funcion <- function(X){
  #1) centrar y reducir X
  X <- centrar_y_reducir(X)
  
  #2) matriz de correlaciones
  X_R <- matriz_correlaciones_ind(X)
  
  #3y 4) valores y vectores propios de la matriz de correlaciones
  X_R_e <- eigen(X_R)
  X_R_valores_propios <- X_R_e$values #ya vienen ordenados de mayor a menor
  X_R_vectores_propios <- X_R_e$vectors
  X_R_vectores_propios
  
  #5) Construir matriz V que tiene como columnas los vectores propios de la matriz de correlaciones
  V <- X_R_vectores_propios

  #6) Calcular la matriz de componentes principales 
  C <- X%*%V 
  
  #7) Calcular la matriz de calidades de los individuos
  X_Q <- matriz_cos2_ind(C, X)
  
  # La matriz de contribuciones para los individuos de X es:
  X_contrib <- contribucion_ind(C, X_R_valores_propios)
  
  #8) La matriz de coordenadas de las variables
  X_T <- coordenadas_var(V, X_R_valores_propios)
  
  #9) Matriz de calidades de las variables
  S <- X_T^2
  
  # Matriz de contribuciones de las variables
  var_contrib <- contribuciones_var(S)
  
  #10) Vector de inercias
  X_I <- calcular_inercias(X_R_valores_propios)
  
  resultado <- list(valores_propios = X_R_valores_propios,
                    ind = list(coord = C, cos2 =   X_Q, contrib = X_contrib),
                    var = list(coord = X_T, cos2 = S, contrib = var_contrib),
                    inercias = X_I)
}

ACP_X <- ACP_funcion(X)
ACP_X$valores_propios

X_ACP <- PCA(X, ncp = 4, graph = FALSE)
X_ACP$eig

ACP_X$ind$coord

X_ACP$ind$coord

ACP_X$ind$cos2

X_ACP$ind$cos2

ACP_X$ind$contrib

X_ACP$ind$contrib

ACP_X$var$coord

X_ACP$var$coord

ACP_X$var$cos2

X_ACP$var$cos2

ACP_X$var$contrib

X_ACP$var$contrib

ACP_X$inercias[1:2]

plot(X_ACP, axes=c(1, 2), choix="ind", col.ind="blue",new.plot=TRUE)

#---------------------------Plano principal individuos--------------------------

plano_principal <- function(matriz){
  #Obtenemos la matriz con las coordenadas de los individuos
  
  ACP_matriz <- ACP_funcion(matriz)
  C <- ACP_matriz$ind$coord
 
  #Convertimos la matriz en un dataframe
  C_data <- as.data.frame(C)
  col_names <- paste("Dim", 1:ncol(C_data))  # Genera nombres como "Dim 1", "Dim 2", etc.
  colnames(C_data) <- col_names

  if(is.null(rownames(C_data))) {
    C_data$individuo <- seq_len(nrow(C_data))
  } else {
    C_data$individuo <- rownames(C_data)
  }

  #Sacamos la inercia
  X_I <- ACP_matriz$inercias
  X_I_redondeados <- round(X_I, 2)
  inercia <- X_I_redondeados[1] + X_I_redondeados[2]
  
  
  #Graficamos
  
  if(nrow(C_data) <= 10){
    individuos <- ggplot(C_data, aes(x = `Dim 1`, y = `Dim 2`)) +
      geom_point(color = "lightblue") +  
      geom_text(aes(label = individuo), vjust = 0, hjust = 0, 
                color = "lightblue") + 
      labs(x = paste("Dim 1 (", X_I_redondeados[1], "%)"), 
         y = paste("Dim 2 (", X_I_redondeados[2], "%)"), 
         subtitle = paste("Inercia = ", inercia, "%"), title = "Plano principal (Individuos)")+
      theme_minimal()
  }else{
    individuos <- ggplot(C_data, aes(x = `Dim 1`, y = `Dim 2`)) +
      geom_point(color = "lightblue") +  
      labs(x = paste("Dim 1 (", X_I_redondeados[1], "%)"), 
         y = paste("Dim 2 (", X_I_redondeados[2], "%)"), 
         subtitle = paste("Inercia = ", inercia, "%"), title = "Plano principal (Individuos)") +
      theme_minimal()
  }
  
  return(individuos)
}

#----------------------------Círculo de correlaciones---------------------------
circulo_correlaciones <- function(matriz){
  #Obtenemos la matriz con las coordenadas de las variables
  
  ACP_matriz <- ACP_funcion(matriz)
  
  X_T <- ACP_matriz$var$coord
  
  #Convertimos la matriz en un dataframe y ajustamos para el gráfico
  X_T_data <- as.data.frame(X_T)
  
  col_names <- paste("Dim", 1:ncol(X_T_data))  # Genera nombres como "Dim 1", "Dim 2", etc.
  colnames(X_T_data) <- col_names
  
  if(is.null(colnames(matriz))) {
    X_T_data$variable <- seq_len(ncol(matriz))
  } else {
    X_T_data$variable <- colnames(matriz)
  }
  X_T_data$`x origen` <- 0
  X_T_data$`y origen` <- 0
  
  #Sacamos la inercia
  X_I <- ACP_matriz$inercias
  X_I_redondeados <- round(X_I, 2)
  inercia <- X_I_redondeados[1] + X_I_redondeados[2]
  
  #Graficamos
  variables0 <- ggplot(X_T_data, aes(x = `x origen`, y = `y origen`)) +
    geom_segment(aes(xend = `Dim 1`, yend = `Dim 2`), 
                 arrow = arrow(length = unit(0.2, "inches")), color = "orange") +
    geom_text(aes(x = `Dim 1`, y = `Dim 2`, label = variable), vjust = -0.5, 
              nudge_y = 0, nudge_x = 0, color ="orange") +
    labs(x = paste("Dim 1 (", X_I_redondeados[1], "%)"), 
         y = paste("Dim 2 (", X_I_redondeados[2], "%)"), 
         subtitle = paste("Inercia = ", inercia, "%")) +
    theme_minimal()
  
  variables <- variables0 + 
    geom_path(data = data.frame(x = cos(seq(0, 2 * pi, length.out = 100)),
                                y = sin(seq(0, 2 * pi, length.out = 100))),
              aes(x, y), color = "black", linewidth = 1, linetype = "dashed") +
    labs(title = "Círculo de correlaciones")
  
  return(list(variables0 = variables0, variables = variables))
}

#---------------------------Grafico Dual----------------------------------------
grafico_dual <- function(graf_ind, graf_var, matriz) {
  # Obtener los graficos y las etiquetas 
  grafico_circulo <- graf_var$variables0
  data_graf_ind <- ggplot_build(graf_ind)$data[[1]]
  data_graf_var <-ggplot_build(grafico_circulo)$data[[1]]
  
  if(is.null(rownames(matriz))) {
    data_graf_ind$label <- seq_len(nrow(matriz))
  } else {
    data_graf_ind$label <- rownames(matriz)
  } 
  
  if (is.null(colnames(matriz))){
    data_graf_var$label <- seq_len(ncol(matriz))
  }else {
    data_graf_var$label <- colnames(matriz)
  }
  
  # Superponer los gráficos
  if(nrow(data_graf_ind) <= 10){
    grafico_final <- grafico_circulo +
      geom_point(data = data_graf_ind, aes(x, y), color = "lightblue") +
      geom_text(data = data_graf_ind, aes(x, y, label = label), 
                vjust = 0, hjust = 0, color = "lightblue")
  }else{
    #grafico_final <- grafico_circulo +
      #geom_point(data = data_graf_ind, aes(x, y), color = "lightblue", size=1)
    grafico_final <- graf_ind + 
      geom_segment(data = data_graf_var, aes(x = x, y = y, xend =  xend, yend = yend), 
                   arrow = arrow(length = unit(0.2, "inches")), color = "orange")+
      geom_text(data = data_graf_var, aes(x = xend, y = yend, label =label), 
                vjust = 0, hjust = 0, color = "orange" )
  }
  
  return(grafico_final)
}

individuos_FM <- plot(X_ACP)

grid.arrange(plano_principal(X), individuos_FM, ncol = 2)

variables_FM <- plot(X_ACP, axes=c(1, 2), choix="var", col.var="orange",
                     new.plot=TRUE)

grid.arrange(circulo_correlaciones(X)[["variables"]], variables_FM, ncol = 2)

dual_FM <- fviz_pca_biplot(X_ACP,col.var = "orange",col.ind = "lightblue")
grid.arrange(grafico_dual(plano_principal(X), circulo_correlaciones(X), X), 
             dual_FM, ncol = 2)

ind_sup_proyeccion <- function(fila, matriz) {
  
 #Obtener medias de las columnas de la matriz
  medias <-c()
  for(i in 1:ncol(matriz)) {
    medias[i] <- mean(matriz[,i])
  }
 
#desviaciones estándar poblacionales de la matriz 
  sd <- c()
  n <- nrow(matriz)
  for(i in 1:ncol(matriz)) {
    sd[i] <- sqrt(((n-1)/n))*sd(matriz[,i])
  }

  #centramos y reducimos la fila y la matriz
  for(i in 1:length(fila)){
    fila[i] <- (fila[i]-medias[[i]])/sd[[i]]
  }
  
  matriz <- centrar_y_reducir(matriz)
  
  #Matriz de correlaciones
  correlaciones <- matriz_correlaciones_ind(matriz)
  
  #Vectores propios
  matriz_e <- eigen(correlaciones)
  V <- matriz_e$vectors
  
  #Se calcula las coordenadas 
  C <- fila%*%V
  
  #Se gráfica esas coordenadas en el plano mediante la conversión de la fila en 
  # un dataframe.
  C_data_ind <- as.data.frame(C)
  col_names <- paste("Dim", 1:ncol(C_data_ind))  # Genera nombres como "Dim 1", "Dim 2", etc.
  colnames(C_data_ind) <- col_names
 
  if(is.null(rownames(fila))){
    C_data_ind$individuo <- nrow(matriz) + 1
  }else {
    C_data_ind$individuo <- rownames(fila)
  }
  
  resultado <- plano_principal(matriz) + geom_point(data = C_data_ind, aes(x = `Dim 1`, y = `Dim 2`),
                                                    color = "red")+
    geom_text(data = C_data_ind, aes(x = `Dim 1`, y = `Dim 2`), 
              vjust = -0.5, hjust = -0.5, color = "red", label = C_data_ind$individuo )
  return(resultado)
  
}

fila <- X[6,]
coordenadasind_sup <- ind_sup_proyeccion(fila, X[-6,])

ACP_indsup <- PCA( X,ind.sup = 6, graph = F)
coordenadasind_sup_FM <-plot.PCA(ACP_indsup, choix = "ind")

grid.arrange(coordenadasind_sup, coordenadasind_sup_FM, ncol = 2)

var_sup_proyeccion <- function(columna, matriz) {
  #se calcula la media y desviación estándar de la columna 
  media <- mean(columna)
  n <- length(columna)
  sd <-sqrt(((n-1)/n))*sd(columna)
  
  #Se obtienen los componentes principales 
  ACP_matriz <- ACP_funcion(matriz)
  C <- ACP_matriz$ind$coord
  
  #Se calculan las correlaciones de la columna con los componentes principales
  coordenada <- c()
  for(i in 1: ncol(matriz)){
    coordenada[i] <-cor(columna, C[,i])
  }
  
  #Se gráfica en el círculo de correlaciones
  circulo <- circulo_correlaciones(matriz)
  graf_circulo <- circulo$variables
  
  #Convertimos la matriz en un dataframe y ajustamos para el gráfico
  coordenada_var <- as.data.frame(t(coordenada))
  col_names <- paste("Dim", 1:ncol(coordenada_var))  # Genera nombres como "Dim 1", "Dim 2", etc.
  colnames(coordenada_var) <- col_names
  coordenada_var$`x origen` <- 0
  coordenada_var$`y origen` <- 0
  
  if(is.null(colnames(columna))){
   coordenada_var$variable <- ncol(matriz) + 1
  }else {
    coordenada_var$variable <- colnames(columna)
  }
  
  resultado <- graf_circulo + 
    geom_segment(data = coordenada_var, aes(x = `x origen`, y = `y origen`, xend = `Dim 1`, yend = `Dim 2`), 
    arrow = arrow(length = unit(0.2, "inches")), color = "red")+
    geom_text(data = coordenada_var, aes(x = `Dim 1`, y = `Dim 2`), 
              vjust = -0.5, hjust = -0.5, color = "red", label = coordenada_var$variable )
  return(resultado)
}

columna <- X[,3]
coordenadasvar_sup <- var_sup_proyeccion(columna, X[,-3])

ACP_varsup <- PCA( X,quanti.sup = 3, graph = F)
coordenadasvar_sup_FM <- plot.PCA(ACP_varsup, choix = "var")

grid.arrange(coordenadasvar_sup, coordenadasvar_sup_FM, ncol = 2)